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Abstract 

i-A I The dynamical mean-field theory (DMFT) is employed to study the correlation- 

Ch ' driven metal-insulator transition in the semi-infinite Hubbard model at half-filling 



and zero temperature. We consider the low-index surfaces of the three-dimensional 
simple-cubic lattice and systematically vary the model parameters at the very sur- 
face, the intra- and inter-layer surface hopping and the surface Coulomb interaction. 
^ . Within the DMFT the self-energy functional is assumed to be local. Therewith, the 

Cn I problem is self-consistently mapped onto a set of coupled effective impurity mod- 

els corresponding to the inequivalent layers parallel to the surface. Assuming that 
V^ . the infiuence of the high-energy Hubbard bands on the low-energy quasi-particle 

^^ I resonance can be neglected at the critical point, a simplified "linearized DMFT" be- 

Qv ' comes possible. The linearized theory, however, is formally equivalent to the Weiss 

molecular-field theory for the semi-infinite Ising model. This implies that quali- 
tatively the rich phenomenology of the Landau description of second-order phase 
transitions at surfaces has a direct analogue for the surface Mott transition. Moti- 
'^ ■ vated by this formal analogy, we work out the predictions of the linearized DMFT in 

C , detail. It is found that under certain circumstances the surface of a Mott insulator 

can be metallic while a Mott-insulating surface of a normal metal is not possible. 
We derive the corresponding phase diagrams, the (mean-field) critical exponents 
and the critical profiles of the quasi-particle weight. The results are confirmed by a 
rS I fully numerical evaluation of the DMFT equations using the exact-diagonalization 

C^ ' (ED) method. By means of the ED approach we especially investigate the non- 

critical parts of the phase diagrams and discuss the U and layer dependence of the 
quasi-particle weight. For strong modifications of the surface model parameters, the 
surface low-energy electronic structure dynamically decouples from the bulk. 



PACS: 71.10.Fd, 71.30.+h, 73.90.+f 



1 Introduction 

The correlation-driven transition from a paramagnetic metal to a paramagnetic Mott- 
Hubbard insulator |jl|, |^ constitutes one of the fundamental problems in solid-state theory. 
The Mott transition is interesting since strong electron correlations lead to low-energy 
electronic properties that are qualitatively different from those predicted by band theory. 

Since it has been recognized that the limit of infinite spatial dimensions {D = oo) 
is a well-defined and meaningful limit also for itinerant-electron models [^ and since the 
invention of dynamical mean- field theory (DMFT) 0, |^ , there is a renewed interest in the 
Mott transition . The DMFT provides an in principle exact description of the transition 
in infinite dimensions. While this is a somewhat artificial limit, the DMFT, as a mean-field 
concept, represents a powerful approach also for the study of finite-dimensional systems. 
Analogous to the Weiss molecular-field theory for localized spin models, the DMFT can 
be expected to give a valuable mean-field picture of the physics of three-dimensional 
itinerant-electron models 0. 

Presumably, the simplest model that includes the essentials of the Mott transition is 
the Hubbard model p|, |^, 0. From the application of the iterative perturbation theory 
(IPT) |ll[| to the D = oo Hubbard model at half-filling and zero temperature, the following 
scenario for the Mott transition has emerged |^: For small Coulomb interaction U, the 
system is a metallic Fermi liquid with a quasi-particle peak at the Fermi energy in the 
one-electron spectrum. As U approaches a critical value Uc2 from below, the quasi-particle 
weight vanishes continuously, similar as in the Brinkman-Rice approach |T^. For strong 



U the system is a Mott-Hubbard insulator with an insulating gap in the one-electron 
spectrum, similar to the Hubbard-III approach |T^. The insulating solution ceases to 



exist when U approaches another critical value Ud from above. In the entire coexistence 
region Ud < U < Uc2 the metallic solution is stable, and thus the transition is of second 
order. The pre-formed gap opens discontinuously at Uc = Uc2- 

It has been questioned [|1^ [15|, |1^, [l^ whether the picture given by the IPT is correct. 



Recent numerical renormalization-group calculations (NRG) |T8|, |T^], however, corrobo- 
rate the IPT results qualitatively, although Uc is found to be significantly smaller than in 
the perturbational approach pO|. On the other hand, there is a remarkable agreement of 



the NRG with the result for Uc in the projective self-consistent method (PSCM) pT| . 

The NRG calculations show that for U ^—>- Uc the quasi-particle resonance becomes 
more or less isolated from the high-energy Hubbard bands |[T^. The resonance basically 
reproduces itself in the self-consistent evaluation of the mean-field equations. This fact 
can be used for a simplified treatment of the DMFT where the influence of the Hubbard 



bands on the low-energy peak is neglected ||2^. This "linearized DMFT" yields a simple 



algebraic equation for the quasi-particle weight at the critical interaction and thereby 
allows for an analytical estimate of Uc- The results are in good agreement with the 
numerical values for Uc obtained from NRG and PSCM on different lattices p2|. For 



inhomogeneous systems, the linearized DMFT also determines the critical profile of the 
quasi-particle weight and the dependence of Uc on the system geometry. Comparing with 
the numerical results obtained from the exact diagonalization method (ED), a convincing 
qualitative agreement with respect to the thickness and geometry dependence of Uc has 
been found for thin Hubbard films 



It has been noticed EH] that the equation of the linearized DMFT is of the same form 



as the linearized mean-field equation of the Weiss molecular-field theory for the Ising 
model (at the critical temperature). There is a one-to-one correspondence if one identifies 
the quasi-particle weight z with the magnetization m, the squared interaction t/^ with 
the temperature T, and the squared hopping integral t^ with the exchange coupling J: 

z^m, f/2 ^ ksT , 36t2 ^ J/2 . (1) 

The Weiss theory, on the other hand, can be considered as being a coarse-grained realiza- 
tion of the classical Landau theory of second-order phase transitions [Q . Consequently, 
the results of Landau theory (for T = T^) can be translated back into predictions con- 
cerning the Mott transition in the Hubbard model (for V = Uc)- 

While the Landau theory of phase transitions is rather simple for homogeneous sys- 
tems, the mean-field theory of critical behavior at surfaces is much more involved and 



numerous non-trivial results can be derived [25|. The idea of the present paper is thus to 



take the Landau theory as starting point and motivation to work out the predictions of 
the linearized DMFT for Hubbard surfaces and finally to test the predictions, as far as 
possible, by comparing with a fully numerical solution of the DMFT equations. 

Within the classical Landau theory, the free energy is expanded in terms of the local 
order parameter m{r). For a semi-infinite system (surface geometry) one additionally 
considers a surface contribution to the free energy |^. Laterally, the order parameter is 
assumed to be homogeneous. We take the x axis be parallel to the surface normal and 
pointing into the volume (x > 0), then m = m{x), and m{x = 0) is the surface value of 
the order parameter. Let us list those mean-field predictions derived from the Ginzburg- 
Landau free energy [^] which - by means of the above-mentioned formal analogy - have 
a direct counterpart for the Mott transition: 

1. The transition in the bulk of the semi- infinite system occurs exactly at the same 
critical temperature Tc^hnik as for the infinitely extended system: Tc^buik = Tc. 

2. Near the surface the order-parameter profile m{x) vanishes at a distance A beyond 
the surface if linearly extrapolated from the boundary. The so-called "extrapolation 
length" A as well as the (bulk) correlation length ^ are the two length scales that 
characterize the order-parameter profile in the continuum model. Microscopically, 
the extrapolation length is related to the model parameters at the surface. In the 
molecular-field approximation of the Ising model we have A~^ oc ( A^ — A) where A 
is the modification of the exchange coupling within the surface layer, Jn = J(l + A), 
and A = Ac corresponds to A = cxd. 

3. For uniform parameters (A = 0) the mean field is smaller at the surface due to 
missing neighbors. This implies a weaker tendency to order. m{x = 0) is smaller 
than m{x i-^ oo) = "mbuik, and m{x) monotonously increases with increasing x (this 
implies A > 0). There is a finite order parameter m(x = 0) > at the surface 
only for T < T^ i. e. only when there is spontaneous order also in the bulk. The 
transition at Tc is termed the "ordinary transition" . 

4. For A < (A > Ac) the surface layer orders at a temperature Tc^surf > ^c,buik ( "sur- 
face transition"). In the regime Tc^buik < T < Tc surf, the bulk correlation length ^ is 
finite, and the order parameter decays exponentially fast from its maximum value 



m{x = 0) at the surface towards zero in the bulk. At T = Tc^buik ("extraordinary 
transition"), the divergence of ^ and the onset of order in the bulk induce singular- 
ities in the behavior of surface response functions. For the order parameter at the 
surface m(x = 0), there is a discontinuity of its second derivative only. At = T^buik 
the order-parameter profile decays algebraically, m{x) oc 1/x. 

5. It holds: (Tc,surf - Tc,huik)/Tc,huik oc A"^. The transition at T = Tc,surf = Tc^hnik in 
the case A = oo is called the "special transition" . For A = oo the order-parameter 
profile is fiat in the ordered phase; the trivial solution m{x) = mtuik = const. 
minimizes the Ginzburg-Landau free energy. In this situation the effect of missing 
neighbors at the surface is compensated exactly. The topology of the phase diagram 
(ordinary, extraordinary, surface and special transition) should be correct whenever 
the surface can support independent order |^. For example, there is no surface 
transition in the semi-infinite two-dimensional Ising model since the "surface" is 



one-dimensional ||27[| . 



6. There are two critical exponents that are merely related to the critical temperatures 
(instead of describing the critical behavior of order parameter and response func- 
tions), the "shift exponent" A^ and the "crossover exponent" 0. They are defined as 
(Tc{d)—Tc{d = oo)) /Tc{d = cxd) oc rf"'^'' for rf i— > oo where Tc{d) is the critical temper- 
ature of a film of finite thickness d, and (Tc^surf (A) — Tc,buik)/7'c,buik oc (A/ Ac —1)^^'^ 
for A 1-^ Ac (special transition). Within Landau mean- field theory one has A^ = 2 
and = 1/2. 

7. Spontaneous order in the bulk always induces a finite order parameter at the surface, 
m{x = 0) > 0. 

The Landau theory also makes additional statements concerning e. g. the bulk and 
surface critical exponents of the order parameter as well as the exponents of response 
functions with respect to an external applied field. We do not mention such results in 
the present context, since either they have no obvious analogue for the Mott transition 
(applied field) or they refer to temperatures T i— > Tc but T ^ T^ where the mean-field 
equation cannot be linearized and where the formal analogy ([I|) breaks down. We will, 
however, discuss a simple extension of the linearized DMFT for U ^-^ Uc but U ^ Uc 
which recovers the result z oc (f/c ~ U) of the PSCM [ pT| . 

To a certain extent, the phase diagram predicted by the Landau theory or, respectively, 
by the linearized DMFT can be tested by comparing with a fully numerical evaluation of 
the DMFT equations. We need an approach that is sufficiently simple for a systematic 
study of a large number of geometries and model parameters. For this purpose the exact- 



diagonalization method of Caffarel and Krauth |^ is well suited. We mainly focus on 
the non-critical parts in the phase diagram where the ED is able to give reliable results. 
Critical exponents, for example, cannot be calculated reliably. The ED has successfully 
been employed beforehand for the discussion of the Mott transition in thin Hubbard films 
2B| and at Hubbard surfaces [^ where the film and surface electronic structure has been 



discussed in detail. Contrary, the present paper focuses on the surface modification of the 
model parameters and on surface phases and thereby substantially extends the previous 
studies. 



The Mott transition at a surface of the semi-infinite Hubbard model has recently been 
investigated in a paper of Hasegawa [^ on the basis of the Kotliar-Ruckenstein slave- 



boson theory [^. With the present study we methodically improve upon Hasegawa's 
work. We will also show that for U ^^ Uc the perturbation of the system that is introduced 
by the presence of the surface deeply extends into the bulk. It is thus insufficient to 
assume (local) physical quantities to be different from their value in the bulk only in the 
first few surface layers. Such a restriction gives rise to unphysical singularities, e. g. in 



the U dependence of the quasi-particle weight as they are seen in Ref. [pn[. Within the 



slave-boson theory it is found that under certain circumstances a metallic surface can 



coexist with an insulating bulk ||30|. Crucial for the existence of this surface phase is a 
considerable decrease of U at the surface. This is an interesting and also plausible result 
although the required strong decrease of U at the surface appears to be quite unrealistic 
for real systems. 

A physically more relevant modification of the model parameters is, in first place, the 
enhancement or decrease of the hopping integrals at the surface. This may be caused by 
a relaxation of the interlayer distance, for example. According to the scaling law t ~ r~^ 
for d electrons (cf. e. g. Ref. |^), a top-layer relaxation Ar/r of a few per cent results 
in a strong change of the hopping integral. A surface modification of t up to about 10%- 
20% appears to be realistic. Besides the hopping we will also discuss a modification of U 
at the surface. In 3d transition metals, however, this effect seems to be less important 



33| , |34| . In any case, U is expected to be enhanced at the surface ||3J]. On the contrary, 
it will be shown that the interesting surface phase occurs for lowered surface U. Another 
important aspect is the surface geometry which is expected to affect the surface phase 
diagram considerably. Open surfaces with a strong reduction of the surface coordination 
number will show the most pronounced surface effects in the electronic structure. We 
thus consider different low-index surfaces oi a D = 3 simple-cubic lattice. 

The basic assumption of DMFT is the strict locality of the self-energy functional. 
For D = 3 dimensions this represents a strong simplification of the problem. The local 
approximation is well justified for the weak-coupling regime, also for the case of surface 
geometries (see the discussion in Refs. p9| , ^ , ^). For the intermediate- to strong- 
coupling regime, however, the assumption may be questioned. One could alternatively 
investigate a surface oi a. D = oo lattice where the DMFT becomes exact. While this 
will be discussed briefiy, we otherwise consider surfaces in D = 3 dimensions. As in Refs. 



23| , ^ we expect the mean- field concept to be a good starting point for D = 3. 



The plan of the paper is the following: The next section introduces into the model. 
The application of DMFT for surface geometries is briefiy discussed in Sec. ^j. We use 
two different methods to solve the DMFT equations: The first one is the approximative 
linearization of the equations for U = Uc p^. This is presented in Sec. ||. Sec. ^ then 
gives a discussion of the analytical results. For the full solution of the DMFT equations, 
we employ the exact-diagonalization method |^ which is introduced in Sec. ^. The 
corresponding results are discussed in Sec. |^. Finally, Sec. |] concludes the paper. 



2 Semi-infinite Hubbard model 

We investigate the Hubbard model on a three-dimensional, simple-cubic and semi-infinite 
lattice. The system is considered to be built up by two-dimensional layers parallel to the 
surface. Accordingly, the position vector to a particular site in the semi-infinite lattice 
is written as Rsitc = i"j + H-a- Here R^ stands for the coordinate origin in the layer a, 
and the layer index runs from a = 1 (topmost surface layer) to infinity (bulk). Fj is the 
position vector with respect to a layer-dependent origin and runs over the sites within the 
layer. In this notation the Hamiltonian reads: 

^ / ^ '^ia,jl3 ^iacr^jPcr i / ^ '''iaa'''ia—a ■ \^ ) 

ijaf3cr icta 

a =t, I is the spin index, f/^ is the (layer-dependent) Hubbard interaction strength. The 
hopping integrals are restricted to be non-zero between nearest neighbors. The energy 
zero is defined by setting tia^ia = to = for sites in the bulk [a. i-^ oo). The energy scale 
is given by taking the (bulk) nearest-neighbor hopping to be t/^iajf^) = —t with t = 1. 

The presence of the surface implies a breakdown of translational symmetry with re- 
spect to the surface normal direction. Lateral translational symmetry, however, may be 
exploited by performing a two-dimensional Fourier transformation: 

ea/3(k) = l^Ee-^'^^'''"^^^-.7^- (3) 

II ij 

Here k is a two-dimensional wave vector of the first surface Brillouin zone (SBZ), and 
A^ll denotes the number of sites within each layer (A^y i-^ oo). Let us briefiy discuss the 
Fourier-transformed hopping matrix which reads: 

/ tiie||(k) + Ato ti2ex(k) N^ 

t2ie±(k) teji(k) tex(k) 

te^(k) te||(k) .. 

V .. .. y 



e(k) 



(4) 



For a > 2 the intra-layer (parallel) hopping and the inter-layer (perpendicular) hopping 
are written as eaaiX) = ^e||(k) and eaa+i(M) = te_i(k), respectively. We consider three 
different low-index surfaces of the sc lattice. The hopping matrix for the sc(lOO) surface 
is obtained from: 

e||(k) = — 2(cos(A;a;) + cos(A;j^)) , 

|ex(k)p = l. (5) 

The perpendicular hopping is k-independent in this case. For the (110) surface we have: 

e||(k) = -2cos(A;a;) , 

\e±{k)\^ = 2 + 2cos{y/2ky), (6) 

and for the sc(lll) surface: 

e||(k)=0, 

|ex(k)p = 3 + 2cos(V2A;j^) 

+ 4 cos{^f3/2k^) cos{^fl/2ky) . (7) 
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Since two nearest neighbors are always located in two different (adjacent) layer, the intra- 
layer hopping vanishes in the last case. Note, that only the absolute square of e_L(k) enters 
the physical quantities we are interested in. 

At the very surface of the semi-infinite system, we consider three different possible 
modifications of the model parameters. Firstly, the intra-layer hopping within the topmost 
surface layer tn may differ from its bulk value (see Eq. (^)). Secondly, we allow for an 
altered hopping ti2 = ^21 7^ ^ between the topmost and the sub-surface layer. Finally, the 
on-site Coulomb interaction strength is assumed to be layer-independent Ua = U = const. 
except for the topmost layer, Ua=i 7^ U . 

We restrict ourselves to the case of manifest particle-hole symmetry, namely a bipartite 
(sc) lattice, nearest-neighbor hopping and half-filling {n = 2{niacr) = !)• In this case the 
Fermi energy is given by /i = to + U/2. It is fixed by the bulk values for the on-site hopping 
and for the Hubbard interaction. Consider the atomic limit t = for a moment: The 
positions of the two Hubbard "bands" in the bulk spectrum are given by: i^iow = to — fi 
and ^'high = to~ f^ + U, i. e. they lie symmetric with respect to /i. In thermal equilibrium 
/x is also the Fermi energy for the top layer. The Hubbard peaks in the surface density 
of states lie at i^iow = to — fJ^ + Ato and -Ehigh = ^0 + ^to — /^ + ^1, where we have taken 
into account the top-layer modification of the interaction strength and where we have 
introduced an additional modification Ato of the atomic level for top-layer sites (see Eq. 
d^)). To maintain manifest particle-hole symmetry and to ensure {jiiaa) = 0.5 for a = 1, 
the Hubbard peaks must again lie symmetric with respect to /z. Thus, we need: 

Ato = (f/-f/i)/2. (8) 

With this choice for the top-layer on-site hopping, the local density of states (DOS) 
Pa{E) = (— I/tt) Im {{ciaa'i Haa)) E IS a Symmetric function of energy for each a. 

We finally introduce the intra- and inter-layer coordination numbers q and p which 
denote the number of nearest neighbors within the same layer and in one of the two 
adjacent layers, respectively. We have: 

g = 4 , p = 1 for sc(lOO) 
g = 2 , p = 2 for sc(llO) 
g = , p = 3 for sc(lll) . (9) 

The bulk coordination number is Z = q + 2p. The surface coordination number is Zs = 

q + p. 

3 Dynamical mean-field theory for surface geome- 
tries 

The one-particle Green function {{ciaa] (^jpu)) contains any important information we are 
interested in. Its diagonal elements Ga{E) = Gia,ia{E) = {{ciaa] cj^a)) e can be written 
in terms of the hopping matrix (H) and the self-energy: 




^°'^'^^^l (^ + .)l-<k)-E(^) ) ■ (^«' 



The self-energy matrix is taken to be k-independent and diagonal, T,ai3{E) = 6ai3^a{E), 
with layer- dependent elements: We assume that the self-energy is a strictly local quantity. 
In the case of an infinitely extended lattice with full translational symmetry, this ba- 
sic assumption leads to the well-known equations of dynamical mean-field theory [^, |^ 
which self-consistently map the bulk lattice problem onto a single-impurity problem 
[|TT], ^. The present case of reduced translational symmetry can be treated analo- 
gously: A local self-energy implies that the Luttinger-Ward functional |^ depends on 
the local (but layer-dependent) propagators only: $ = $[..., Ga{E), . . .]. This in turn 
means that the self-energy of the a-th layer is solely a functional of the local propagator: 
Ec(£') = S^/SGa{E) = S[Ga{E)]. The functional S is the same as in the case of an im- 
purity problem, e. g. the single-impurity Anderson model (SIAM), Simp(-E) = S[Gimp{E)], 
because the same type of skeleton diagrams occur in the expansion of the impurity self- 
energy Simp(-E'). With each layer a = 1, 2, ... we therefore associate a SIAM, 

+ E W"^4.c. + H.c.) , (11) 

ka 

with e|f = Ua^ia and where the conduction-band energies e)^ and hybridization strengths 
V^ chosen such that we have 

A(")(E) =E + f,-e'f^- J:1:;^{E) - G^{E)-' (12) 

for the hybridization function A^"\E—h) = J2ki^k y/i-^~^k )■ i^^- (HD ^^^Y provides 
an implicit definition of the hybridization function since S^^ depends on A*^")). This im- 



(") 



plies at once the equahty between the impurity Green function of the a-th SIAM, G\^p{E), 
and the on-site lattice Green function in the a-th layer Ga{E) and thus the equality be- 
tween the respective self-energies, S>^p(i?) = S[G\^p{E)] and T^aiE) = S[Ga{E)]. 

The following iterative procedure then allows to solve the semi-infinite Hubbard model 
within the dynamical mean-field approximation: Starting from a guess for the layer- 
dependent self-energies Tja{E), we calculate the on-site Green function of the a-th layer 
using Eq. (p^Of). Via Eq. (0), Ga{E) and Sa(-E') = T.l'^Jp{E) determine the hybridization 
function A^^'^E) of the a-th SIAM. The crucial step consists in solving the impurity 
models for a = 1,2,... to get the impurity self-energies S>°p(ii^) which are required for 
the next cycle. The cycles have to be repeated until self-consistency is achieved. 

Applying the DMFT to the semi-infinite Hubbard model means to map the original 
lattice problem onto an infinite set of impurity problems. The mapping is mediated by 
the self-consistency equation ([12|) for a = 1,2, ...,cx). For a given set of hybridization 
functions, each impurity model can be treated separately. There is, however, an indirect 
coupling via Eq. (|I0|) in the self-consistency cycle since the on-site Green function of 
a given layer depends on all layer-dependent self-energies. The essential difference with 
respect to the case of an infinitely extended lattice with full translational symmetry where 
only one single-impurity model and only one self-consistency condition is needed, consists 
in the fact that for a semi-infinite system the sites within different layers have to be 
considered as non-equivalent. 
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4 Linearized DMFT at the critical interaction 

The zero-temperature Mott transition from a paramagnetic metal to a paramagnetic 
insulator is actually hidden due to antiferromagnetic order which is realized in the true 
ground state. To study the Mott transition, the solutions of the mean-field equations have 
to be enforced to be spin symmetric. There have been numerous DMFT studies of the 
T = Mott transition in the recent past using different methods to solve the impurity 



problem: the iterative perturbation theory (IPT) [0, ^ |^, the exact-diagonalization 
approach (ED) |2^, ^ ^, the projective self-consistent method ||2T| as well as numerical 
renormalization-group calculations (NRG) |T8|, |19 . 

The IPT and, in first place, the NRG results show that for ?7 i— > f/^ the quasi-particle 
resonance at ii^ = is more or less isolated from the high-energy Hubbard peaks at 
E ^ ±f//2. The resonance basically reproduces itself in the self-consistent procedure to 
solve the DMFT equations. A reasonable assumption is therefore that for U = Uc the 
low-energy part of the SIAM hybridization function A^"\E) consists of a single pole at 
E = only, 

a(") 

A(")(i^)^-^, (13) 

and that the effect of the Hubbard bands can be disregarded completely. With this 
assumption a simplified, "linearized" DMFT becomes possible P^ , p5| . There is an at- 
tractive feature of this method which outweights the necessity for a further approximation: 
It allows for a fully analytic treatment of the mean-field equations, and an analytical ex- 
pression for Uc is obtained. Studying the dependencies of Uc on the model parameters 
can provide a valuable first insight into the problem. The predictions of the linearized 
theory have been compared beforehand with fully numerical DMFT results for the Bethe 



and the hypercubic lattice in D = cxd [^ and for the case of thin Hubbard films ||23 
A satisfactory quantitative agreement has been noticed. This makes us confident that at 
least the correct trends can also be predicted for the case of a semi-infinite lattice. 

The details of the method can be found in Ref. |2^; here we simply repeat the main 
idea and the final result: In the ansatz for the hybridization function (|T3|), A^ denotes 
the layer-dependent coefficient in the A^-th step of the self-consistency cycle. The aim 
is to calculate A]^|^. The one-pole structure of the hybridization function corresponds 
to a well-defined SIAM with n^ = 2 sites which can analytically be solved for each a. 
In the one-particle excitation spectrum of the a-th SIAM there are two 6-peaks a.t E ^ 
±Ua/2 as well as two (5-peaks near E = corresponding to the (infinitely sharp) Kondo 
resonance for U = Uc in the infinite (n^ = cxd) system. The layer-dependent weight 
of the resonance Za can be read off from the solution. Za determines the self-energy 
T,a{E) = Ua/^ -|- (1 — z~^)E + ■ ■ ■ and via Eq. (p!0D the on-site Green function of the 
a-th layer at low energies. Using these results in the self-consistency equation (|T2]) and 
insisting on the one-pole structure of the hybridization function, yields a new coefficient 
A}^^_j^. At this point the possible infiuence of the Hubbard bands is ignored. The final 
equation that relates Aj^l^i to /Sr^' reads: 



where we have defined the following semi-infinite tridiagonal matrix: 

( qt\jUl ptut/UiU \ 

ptut/UiU qtyu^ ptyw^ 

ptyu^ qtyu^ .. 
\ .. .. J 



K = 36 



(15) 



A self-consistent solution of the linearized mean-field equation (p!4D is given by a non- 
trivial fixed point of K. Let A,- denote the eigenvalues of K. We can distinguish between 
two cases: If |Ar| < 1 for all r, there is the trivial solution limN^oo '^n = only. This 
situation corresponds to the insulating solution beyond the critical point. Contrary, if 
there is at least one A,. > 1, A]^ diverges exponentially as A^ ^— i> oo. This indicates the 
breakdown of the one-pole model for the hybridization function in the metallic solution 
below the critical point. The maximum eigenvalue thus determines via 

Amax = Amax(g, P, U, tn, tu, Ui) = 1 (16) 

the critical model parameters. 

At the critical point the mean-field equation (0) can be written as Za = J^/s KapZp 
since the layer-dependent quasi-particle weight z^ oc A^"). Formally, this equation can 
be compared with the Weiss mean-field equation for the layer magnetizations nia in the 
semi-infinite Ising model with coupling constant J . The linearized mean-field equation for 
T = Tc reads: m^ = {J/2kBT){qma +pma+i +P^a-i) (we assume the model parameters 
at the surface to be unmodified for the moment). The formal analogy with Eqs. (|T^ 
and (|15D is obvious and justifies the identification made in Eq. (P and the corresponding 
discussion in Sec. |l|. 

5 Analytical results 

From the basic equation (^) we can calculate the critical parameters for different cases. 
First, we consider a system that is built up by a finite number of d layers (film geometry). 
The model parameters are taken to be uniform, i. e. tn = ti2 = t and Ui = U (at both 
surfaces) . The eigenvalues of the (i-dimensional matrix (|T3|) can be calculated analytically 



for this case [Q: 



with r = 1, ...,d. Taking the maximum eigenvalue and solving for U yields the thickness 
dependence of the critical interaction ^^ : 



U,{d) = QtJq + 2pcos(^^^ 

Expanding the result for Uc in the limit oi d \-^ oo yields: 

U,{d) - Ucjoo 

with a "shift exponent" A^ = 2. 
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oc d-^' (19) 



In the limit d = oo any of the two film surfaces represents a semi-infinite system. From 
T8| ) we get for the critical interaction: 



f/c,buik = mJq + 2p = 6ty/Z, (20) 



which is the same result as is found when applying the method to the infinitely extended 
bulk system directly |2^. We notice that for the case of uniform model parameters the 
linearized DMFT yields a unique critical interaction for the semi-infinite system which is 
the same as the bulk value. No surface phase is found. This observation is fully consistent 
with what has been obtained in previous numerical DMFT studies of the Mott transition 
at Hubbard surfaces [^ for uniform parameters. Despite the fact that at the surface 



the electronic structure has turned out to be modified considerably, a surface critical 
interaction different from the bulk value has not been found. 

In the following we thus concentrate on a semi-infinite system with modified parameteTs 
at the very surface. Also in this case the condition ([l6|) can be treated analytically: To 
simplify the notation let us write Ku = a', K12 = K21 = b' and Kaa = a, Kaa±i = b for 
a > 2 {a,b,a',b' > 0). Let K(n) be the matrix that is obtained from the semi-infinite 
matrix K = K(0) by deleting its first n rows and columns. Furthermore, we define 
g^{X) = det(Al - K{n + l))/det(Al - K(n)). ^„(A) is the (1,1) or "surface" element of 
the Green matrix (Al — K(n))~^. Expanding the determinant det(Al — K(n -|- 1)) with 
respect to the upper left element, one easily verifies the recurrence relation Qn{^)~^ = 
A — a — b'^Qn+iW for n > 1. However, all the ^n(A) for n > 1 must be equal since the 
(off-)diagonal elements of K(?2 > 1) are constant. This results in a quadratic equation for 
Q, the solution of which is given by: 

^(A) = ^„>i(A) = ±i^X-aT ^(A - a^ - 4b^^ (21) 

for ±(A — a) > 0. The eigenvalue spectrum of the semi-infinite matrix K consists of a 
continuous bulk part which can be read off from Eq. (|T^ for (i t-^ 00 to be given by: 

|A-a|<26. (22) 

This is just the region where Im^(A) 7^ 0. The largest eigenvalue in the bulk continuum 



is given by A = a + 26 corresponding to the bulk critical interaction given in (pO]) . At 
U = f/c,buik the bulk undergoes the metal-insulator transition irrespective of the state of 
the surface. 

Under certain circumstances an additional discrete ( "surface" ) eigenvalue A^ may split 
off the bulk continuum. If a discrete eigenvalue exists, we must have Qo{^s)~^ = 0. Using 
the result (^1]) to determine Qo{X) from the recurrence relation GoW^^ = \ — a' — b''^Q{X), 
we obtain the following equation for the eigenvalue: 

A. - a' - b" (^^ (a, - a T \/(A,-a)2-462)) = . (23) 

(±(As — a) > 0). Solving the equation for A^ yields the position of the eigenvalue in the 
spectrum of K. Since K is real and symmetric, only a solution A^ with ImA^ = is 
meaningful; a discrete A^ must lie outside the bulk continuum (^). Only the maximum 
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eigenvalue in the spectrum is physically relevant [Eq. ([T^)]. Thus, we are interested in a 
solution that is split off the upper edge of the continuum: 

Re A, > a + 26 . (24) 

Since 6 > only the — sign must be considered in (p^). 

Whether or not the condition (p^ can be met, depends on the (surface) parameters 
a' and h' . Solving Eq. (^) for As and inserting the solution into (P^, yields the following 
relation for a' and h' 

262 ^ ^(^ _ ^/) _ 5'2 ^ Q ^ ^25) 

which must be fulfilled to obtain a (physically relevant) surface mode. Note that the 
relation cannot be satisfied with uniform parameters, i. e. a' = a and 6' = 6. 

The interpretation is the following: In a semi-infinite system with surface parameters 
til, ti2 and f/i that do not obey the condition (P3|), there is only the "ordinary" transition 
from a metallic to a Mott insulating state atU = f/c,buik when increasing the interaction 
strength. The critical interaction [/c,buik is given by Eq. (pOf). At this point all layer- 
dependent quasi-particle weights Zq,, in the bulk as well as at the surface, vanish. On 
the other hand, for a sufficiently strong modification of tn, ti2 or f/i, i. e. for a' and 
6' satisfying (pS]), there are two critical interaction strengths: The first one is f/c,buik 
again. At f/ = f/c,buik the bulk quasi-particle weight ^0=00 vanishes. The second critical 
interaction strength t/c,surf can be determined from A^ = 1 where A^ is the solution of 
([23|). Let us assume that ?7c,surf > f4,buik- For U > [/c,surf the entire system is in the Mott 
insulating phase. For f/c,buik < U < t/c^surf, however, the bulk is a Mott insulator while 
the surface is still metallic. We call the transition sXU = ?7c,buik the "extraordinary" and 
the transition at t/ = f/c,surf the "surface transition" in analogy with the terminology for 



magnetic phase transitions at surfaces p5 



The remaining question is whether or not t/c,surf < t4,buik can be possible. In such 
a situation we would have a quasi two-dimensional Mott insulator on top of a metallic 
bulk for interactions ?7c,surf < U < f/c,buik- However, this possibility is ruled out: Eq. (p^ 
can be rewritten as A^ > t^cbuik/^^- Furthermore, at the critical point U = f/c,surf the 
value As = 1 fulfills Eq. (pSj). But this implies 1 > f^c bulk /^c surf- ^^ ^^^ state that the 
linearized theory predicts that a metallic surface coexisting with a Mott insulating bulk 
is possible while the opposite scenario cannot be realized. 

Arguing physically, if (at the Fermi edge) there is a finite (local) density of states in 
the second and all subsequent layers, this must always induce a non-zero, though possibly 
low density of states in the top layer, and thus an insulating surface phase is excluded: 
Consider the free-standing two-dimensional layer at an interaction strength Ui being suf- 
ficiently strong to force the system to the insulating phase. Let the monolayer be coupled 
to the second and all subsequent layers by switching on the hopping between the top and 
the second layer ti2 7^ 0. If ti2 is finite but too small, the low-energy bulk excitations 
cannot propagate into the top layer and are refiected at the Hubbard gap. However, 
virtual hopping processes are possible which cause (an exponentially damped) weight of 
bulk excitations in the top layer. The exponential damping becomes unimportant in this 
case since it is effective in one layer only. 

For the opposite case of a metallic surface on top of a Mott insulator, however, it does 
become essential: Low-energy excitations can propagate within the surface region since 
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U < Uc,snrf- Because U > ?7c,buik, they cannot propagate into the bulk but are reflected at 
the (bulk) Hubbard gap. While virtual processes always generate some non-zero spectral 
weight at the Fermi edge in each layer, the weight is infinitesimally small asymptotically, 
for a I— i> oo. 

Since critical fluctuations spread out all over the system at a second-order critical 
point, different parts of a system should undergo the transition at a common and unique 
critical value of the external control parameter. The exponential damping of low-energy 
excitations over large distances explains why there can be two critical interactions. This 
is analogous to the case of magnetic phase transitions at surfaces: In a system where a 
magnetic surface coexists with a paramagnetic bulk, the layer magnetization must decay 
exponentially when passing from the surface to the crystal volume. Vice versa: A magnetic 
bulk always induces a finite magnetization in the top layer. The exception is the somehow 
artificial case where the top layer is completely decoupled from the rest system (e. g. 
tl2 = 0). 

5.1 Modified intra- layer surface hopping 

Some more aspects of the metallic surface phase shall be addressed in the following. In 
particular, to discuss the effects of the surface geometry, we refer to the different low- 
index surfaces of the sc lattice mentioned above. Furthermore, it is helpful to consider 
the different types of surface modifications separately. 

We start by considering a modified intra-layer hopping in the top layer: tu ^ t. We 
have: 

a'=^^a, b' = b. (26) 

From Eq. (|25| ) we can deduce that there are two critical interactions provided that 

(27) 




According to ([16|) and (^), the critical interaction strength at which the surface transition 
takes place, is given by: 



Urfiurf — 6t 



a2 ^2 



\lf + -72-^2- (28) 

\ t^ q tf 1 - t^ 



The corresponding phase diagram is shown in Fig. |l]. 

For the (111) surface there is no intra-layer hopping at all (g = 0). A rather moderate 
enhancement of tu (about 12%) is sufficient to obtain a metallic surface phase for the 
sc(lOO) surface. In the case of the sc(llO) surface a stronger modification is necessary. 
These trends are plausible: Obviously, for both surfaces a larger tu means that electrons 
in the top layer are more itinerant and thus tends to delay the transition to the insulating 
state as U is increased. A smaller intra-layer coordination number q counteracts this 
mechanism. Consequently, one needs a stronger enhancement of tu for the (110) surface. 
The ?7-range where a metallic surface coexists with an insulating bulk quickly increases 
at til is increased. For tu \-^ oo one would expect that the energy scales relevant for the 
bulk become meaningless and that the electronic structure of the top layer decouples from 
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Fig. |T]: t\\-U phase diagram as obtained from the hnearized DMFT. 
For \J < C/c,buik the system is metalhc. For U > C/c,buik the bulk is a 
Mott insulator, the surface can be either insulating (left to the phase 
boundary) or metallic (right). Phase boundaries for the (100) and (110) 
surface of the sc lattice. Energy units: nearest-neighbor hopping t = 1. 
Free band width W = 12. 



the rest system. This is predicted correctly by Eq. ( p8|) which yields f/c,surf = Stuy^ in 
this limit, i. e. the critical interaction strength of a free-standing two-dimensional layer. 



5.2 Modified inter-layer surface hopping 

For a modified inter-layer hopping between the top and the sub-surface layer ti2 ^ t we 
have: 

(29) 



a' = a, h' = ^-h. 
A metallic surface of a Mott insulating bulk is possible for 

tl2 > V^t, 



(30) 



irrespective of the type of the surface. The critical interaction strength for the surface 
transition is given by: 



U, 



c,surf 



6t- 



\ 



q + p- 



t^ 
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t' - lt\, - t^ 



(31) 



Fig. ^ shows ?7c,surf as a function of ti2 for the different surfaces. 
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Fig. |2|: t\2-U phase diagram. 



An enhancement of ti2 again means to enhance the itinerancy of electrons at the 
surface. Hopping processes between the topmost and the sub-surface layer become more 
likely. A modification of about 19% is sufficient to suppress the transition to the Mott 
insulating phase at the surface for U > f/c,buik- The surface critical interaction strength 
C^.surf up to which the metallic surface phase persists for a given ti2, is the largest for the 
sc(lll) surface since here the perpendicular hopping is favored by the comparatively high 
inter-layer coordination number p = 3 anyway. In the limit ti2 ^—>- cxd the first two layers 
of the surface will decouple from the bulk. The surface critical interaction strength in 
this limit should be the same as for a bi-layer system with strongly anisotropic hopping. 
Consider, for simplicity, the sc(lll) surface where g = 0. In this case all sites in the 
bi-layer system have the same coordination number p and the bulk formula (^) may be 
applied accordingly. This yields: ^c,surf = Qti2\^, which is consistent with the ti2 ^-* cxd 
limit of Eq. {M). 



5.3 Modified surface Coulomb interaction 

Finally, we consider a modified Coulomb interaction in the top layer, Ui ^ U . In this 
case: 



a 



h' 



-b. 



(32) 



As in the two other cases, we could fix the surface model parameters, vary U and ask for 
the critical interaction strength f/c,surf- For the present case, however, it appears to be 
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more intuitive to consider the bulk f/ to be a fixed quantity and to vary Ui. 

For U above the bulk critical interaction ?7c,buik the bulk is a Mott insulator. The 
system may then become critical with respect to Ui, provided that 



]j q + 2p 
The surface transition takes place at Ui = Uic surf with 



(33) 



f/i 



c,surf 



\ 



f/2 + 36gt2 



\ 



'f/2-36gt2 



- 36V^^ 



(34) 



for U > f/c,buik- Fig. 1^ shows the corresponding phase diagram. For f/ i-^ oo we get 
^i,c,surf = ^t\/(i- This is the critical interaction strength of the free-standing monolayer. 

The results for modified surface Coulomb interaction can be compared with Hasegawa's 
slave-boson approach [^. Qualitatively, the respective U-Ui phase diagrams for the 
sc(lOO) surface look similar. The critical interactions predicted by the slave-boson method 
are somewhat larger compared with the DMFT results. This is typical for the slave-boson 
method [§|. An important difference is found with respect to the "special transition" at 
the tri-critical point U 



f4,buik, Ui = Ui^c = Jiq + p)/iq + '2p) Uc,huik- The linearized 
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Fig. H: U-Ui phase diagram. 
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DMFT predicts 

^IjCjSurf f-^Cjbulk 



u, 



oc 



c,bulk 



for Ui H-^ Ui^c with a "crossover exponent" ( 
found for modified surface hopping tu or ti2: 

^c,surf L^c,bulk / 
OC 



U, 



c,bulk 





1 




(35) 


!) = 1/2. 


The same 


crossover 


exponent is 


^11(2) 
^ll(2),c 


1 




(36) 



where tn^c and ti2,c are defined by the r.h.s. of Eqs. (|27|) and (^), respectively. This 
follows from a direct calculation and can also be seen in Figs. |l], ^ and |^. Contrary, 



within the slave-boson theory of Ref. |30[, f/i,c,surf seems to be constant as a function of 



[/, and a crossover exponent cannot be defined. 

5.4 Profiles of the quasi- particle weight 

The mean-field equation of the linearized DMFT, A]^^_^ = X]/3 -^a/3 ^kr , has a non-trivial 
solution only at a critical point for the Mott transition, e. g. at f/ = ?7c,buik or f/ = t/c,surf- 
This solution is a fixed point of the matrix K, A^^ = limTVi-^oo A}^ , and can be calculated 
as the eigenvector of K belonging to the eigenvalue A = 1 [Eq. (|T6|)]. Since Za oc /SS°'\ the 
eigenvector has the meaning of the critical profile of the quasi-particle weight, i. e. the 
a-dependence of z^ in the limit Zq, i-^ 0. It is uniquely determined up to a normalization 
constant. 

The upper left part of Fig. | shows the critical profile at the sc(lOO) surface for different 
values of tu and U = Uc,hnik or U = t/c,surf, respectively. Za has been normalized to its 
top-layer value zi. For unmodified surface hopping tu = t, the profile is linear. In fact, 
the ansatz z^ oc a solves the mean-field equation z^ = {36t'^/U'^){qZa +pza+i +pza-i) for 
U = Uc^huik = Qty^q + 2p. Physically, this means that at the critical interaction the surface 
effects extend into the bulk up to arbitrarily large distances. Note that this implies that 
actually an infinite number of inequivalent surface layers has to be considered in a fully 
numerical evaluation of the DMFT. 

For U close to f/c,buik but U < f/c,buik, one would expect that the profile converges to 
a finite bulk value: lima^oo ^a = 2;buik > 0. In its present form, however, the linearized 
DMFT is not applicable here. One may consider the following extension |^ of the 



mean-field equation (for simplicity, we discuss the case of uniform model parameters, the 
generalization for modified surface hopping ot Ui ^ U is straightforward): 

Za = -JJY {qZa + pZa+1 + PZa^l) - C zl . (37) 

A quadratic term in Za with a constant coefficient c > has been added. The constant c 
can be fixed by the value for Zbuik or for zi (Ref. [^ yields the explicit value c = 11/9 but 
we do not need the result here). This extension of the linearized DMFT is in the spirit 
of Landau theory, we simply consider the next term in an expansion with respect to the 
"order parameter" Za- As in the Landau theory, higher-order terms in Za or quadratic 
terms that couple the different layers are still neglected. The additional term in Eq. (pT]) 
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Fig. ^: Profiles of the quasi-particle weight for the sc(lOO) surface. Upper left: Profiles 
for different tu. U = C/c,buik for in < 1.25 and U = C/c,surf for *ii ^ 1.25. Upper 
right: tu = t = 1 and different U close to C/c,buik- Lower left: tu = 1.5 and different 
f7/t/c,buik- Lower right: f//C/c,buik = 0.99 and different tf2- The profiles are normalized to 
the top-layer value (r.h.s.) or the bulk value (l.h.s.), respectively. 



ensures a linear t/-dependence of the quasi-particle weight in each layer: Za oc {Uc — U) 
for U I— >^ Uc- This is consistent with the (bulk) critical behavior found within the PSCM 

Using Eq. (pTl) we have calculated the profile of the quasi-particle weight for tu = 1 
and different U < Uc,hu\k, see Fig. ^ (upper right). For t//[/c,buik = 0.9 the quasi-particle 
weight significantly differs from the bulk value in the first few layers from the surface only. 
As U *—>■ t/c,buik, however, the linear trend of Za clearly develops. 

A linear trend of the critical profile is also observed for slightly enhanced surface 
hopping, t1^ = 1.1 and tf-^ = 1.2 (Fig. ||, upper left). For a surface hopping tu = 
Jl +p/q = J 5/ 4: we get the so-called special transition (cf. Eq. (pTl) and Fig. 0). At the 
critical interaction the profile is a constant (Fig. ^, upper left). In this case the effect of 
missing neighbors at the surface is exactly compensated by the enhancement of tu- 

For til > v 5/4 there are two critical interactions, ?7c,buik and Uc^smi- For U = f/csurf, 
Za/zi is at its maximum in the top layer and exponentially decays as a i— >^ oo (Fig. 0, 



upper left). For U < f/c,surf [according to Eq. (^7\)] the decay becomes slower until the 
profile converges to a finite bulk value for U < Uc^hnik (lower left). 

Finally, the lower right part of Fig. ^ shows the profile of the quasi-particle weight 
obtained from Eq. (^) for f//f/c,buik = 0.99 and modified inter-layer surface hopping ti2- 
For t^2 < V^ the profile is a monotonously increasing function when passing from the 
surface to the bulk. t^2 = v^ marks the special transition [see Eq. (|30|)]. Here the profile 
would be constant for a > 2 and U = f/c,buik as can be seen from the mean-field equation 
of the linearized DMFT. For t^g > V^ the quasi-particle weight is enhanced at the surface 
and monotonously decreases for a > 2. 

5.5 Infinite dimensions 

Dynamical mean-field theory rests on the local approximation for the self-energy func- 
tional. Since it is known that in the limit of high spatial dimensions Z^ i-^ oo [^ the 
local approximation becomes exact [H|, it may be interesting to discuss the (somewhat 



artificial) case of a surface of the infinite-dimensional hyper-cubic lattice. 

A Z)-dimensional hyper-cubic lattice may be thought to be built up from {D — 1)- 
dimensional "layers" perpendicular to a D-dimensional spatial direction characterized 
by the set of Miller indices [xi,X2, . . . ,xd]- Cutting the hopping between two adjacent 
layers, one obtains a "(xi,X2, . . . ,X£)) surface". Consider the low-index directions with 
Xi = ■ ■ ■ = Xr = 1 and Xr+i = ■ ■ ■ = xd = 0- For a given site there are q = 2D — 2r nearest 
neighbors within the same layer and p = r nearest neighbors in each of the adjacent layers 
{Z = q + 2p = 2D). 

For r = 1, i. e. a (1000...) surface, a site in the topmost layer has Zg = q + p = 2D — 1 
nearest neighbors to be compared with Z = 2D in the bulk. For D h^ oo the local 
environment of the surface sites is essentially the same as in the bulk, surface effects 
become meaningless. With the usual scaling of the hopping t = t*/\^2D [^, the free 
top-layer local density of states (DOS) is a Gaussian p(°)(E) = exp[—{E/t*Y/2]/{y/2nt*) 
- as in the bulk. 

For r = D one obtains the open (1111...) surface. The surface coordination number 
is reduced to Zs = p = D. This implies a ratio As/A = Zs/Z = (q + p)/{q + 2p) = 0.5 
between the variances of the top-layer and bulk DOS. The results of a simple numerical cal- 
culation are shown in Fig. |^. We notice a strongly modified and strongly layer-dependent 
DOS near the surface which slowly converges to the bulk Gaussian DOS for a ^-^ oo. In 
many respects the results resemble the DOS at the D = 3 sc(lll) surface, in particular 
the oscillation of p^"-' {E = 0) as a function of a [^ . 



In infinite dimensions dynamical mean-field theory is exact also for the semi-infinite 
model. The scaling of the hopping implies G\j ~ 1/vD for the free propagator between 
arbitrary nearest-neighbor sites i and j, and the proof that the self-energy is local, is 
essentially unchanged (see Refs. [0, |^, 0)- The simple linearized DMFT can be developed 



as in Sec. ^. We only have to insert the general expressions for the coordination numbers 
q = 2D — 2r and p = D, and to perform the limit D i— »• oo in the equations ( P7D to (^), 
paying attention to the scaling of the hopping. 

Varying r we can then pass continuously from the most closed (r = 1) to the most 
open (r = D) surface geometry. Consider, for example, a modified intra-layer surface 



hopping. A surface phase is predicted to be existing for tj^ > t*Jl + r/{2D — 2r) [cf. 
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Fig. |: ?7 = layer-dependent density of states pa (E) at the (1111...) 
surface of the D = oo hyper-cubic lattice. Scaled hopping t = t* jyJlD 
with t* = 1. "1" stands for the topmost surface layer, "2" denotes the 
sub-surface layer, etc. 



Eq. (|27|)]; i. e. for all t\^ > t* in the case of the closed r = 1 surface and not at all for 
the r = D surface. For r = 1 the surface critical interaction is given by f/c,surf = Qtli 
[Eq. (pSDJ to be compared with the bulk critical interaction Uc^hnik = 6t* [Eq. (pO])]. With 
increasing r, f/c,surf decreases until f/c,surf = f^c,buik for r = D. 

The other cases may be discussed accordingly. Upon taking the limit D i-^ oo, we 
always obtain non-trivial and plausible results. The discussion is analogous to the D = 3 
case. We conclude that the semi-infinite Hubbard model remains non-trivial for D = oo 
and provides a useful framework for investigating the surface phase. In principle, this 
can be done without approximations by employing the DMFT. Recall, however, that the 
linearized DMFT is still approximate (Sec. ^ and Ref. ||22|| ) . 



6 Exact-diagonalization method 

For a complete numerical solution of the mean-field equations at finite temperatures one 
may employ the Quantum-Monte-Carlo method [B^, H^, H^ 



For T 



the Exact- 



Diagonalization (ED) approach p8|, ^ ^ can be applied and is chosen here. The main 
idea is to map onto a SIAM with a finite number of sites n^. Lanczos technique [Q is used 
to calculate the ground state as well as the T = impurity Green function and self-energy. 
The DMFT equations are solved on the discrete mesh of Matsubara energies where the 
(fictitious) inverse temperature (3 introduces a low-energy cutoff. Details of the method 
can be found in Ref. P|. The surface geometry can be simulated by a slab consisting of 
a finite but sufficiently large number of layers d (for U ^ Uc) ■ The numerical effort then 
increases linearly with d at least. In Refs. |2^, ^ we have discussed the application of 
ED to film and surface geometries. 

ED is able to yield the essentially exact solution of the mean-field equations in a 
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parameter range where the errors introduced by the finite system size are unimportant. 
For the Mott problem the relevant low-energy scale is set by the width of the quasi-particle 
peak in the metallic solution. It has to be expected that there are non-negligible finite- 
size effects when this energy scale becomes comparatively small. We are thus limited to 
interactions strengths that are not too close to t/c,buik or ^c,surf and cannot access the very 
critical regimes. This also implies that a precise determination of t/c,buik and f/c,surf and 
thereby a direct comparison with the linearized DMFT is not possible. The discussion in 
[|23|| , however, shows that the main trends can be derived safely. 

In the following we mainly focus on the low-energy electronic structure which the 
ED method is able to predict reliably in the non-critical regimes. The so-called layer- 
dependent quasi-particle weight, 



'^ dE^jE = 0) 



dE J ' (^^) 

is the primary quantity of interest. 2;^ < 1 is weight of the coherent quasi-particle peak 
in the local DOS Pa{E) of the a-th layer or, alternatively, the reduction factor of the 
discontinuous drops in a-th momentum-distribution function ^^(k) when k crosses the 



one-dimensional Fermi "surfaces" 23 



Routinely, the calculations have been performed with Ug = 8 sites in the effective 
impurity problems. For the fictitious temperature we have chosen /3~^ = 0.0016 W {W = 
12 is the free band width). Ug and f3 determine the "energy resolution" which is found 
to be about AE = 0.12 = VT/IOO. This implies that reliable results can be expected in 
a parameter region where Za > 0.01 (cf. Ref. |2^). A moderate number c? < 25 of layers 
in the slab is sufficient to simulate the semi-infinite system - except for the very critical 
regime. This has been checked by comparing the results from calculations for different 
d. We made use of the mirror symmetry at the center of the slab and of electron-hole 
symmetry to reduce the number of parameters, the conduction-band energies ek and the 
hybridization strengths Vk [k = 2, ...,71^), which have to be determined self-consistently. 
We always found a unique and fully stabilized solution. 

7 Numerical results for the sc(lOO) surface 

To keep the calculations manageable, we restrict the discussion to the D = 3 sc(lOO) 
surface in the following. We start with the case of uniform model parameters. Fig. ^ 
shows the bulk quasi-particle weight z (dashed line) as a function of U. It starts from 
its non-interacting value z = 1. A quadratic U dependence is noticed for small U in 
agreement with perturbation theory p9|. z vanishes as U approaches f/c,buik- The overall 



dependence on U is very similar to what is known from DMFT studies of the D = 00 
Bethe lattice §. 

In the top layer of the sc(lOO) surface the quasi-particle weight is significantly reduced 
(solid line). The lowered coordination number at the surface implies a reduced variance As 
of the free surface DOS and thereby an increased effective interaction U/^/As compared 
with the bulk. Thus, at the surface correlation effects are enhanced, and Za=i is lowered. 
Despite this tendency towards an insulating surface, we find a common critical interaction 
Uc,snrf = ^c.buik which, for uuiform parameters, is in agreement with the analytical results. 
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f^c,buik also represents the critical interaction for all sub-surface layers. For the rather 
closed (100) surface, Za{U) is almost identical with the bulk function for a > 2. 

From Fig. | we can read off t/c,buik ~ 16.0 while Eq. (|^) predicts t/c,buik = 14.7. We 
have to bear in mind, however, the underlying assumptions that lead to (^). Moreover, 
as concerns the ED, finite-size effects prevent a precise estimate: ^c,buik ~ 15.1 is found 
for Hg = 10 sites in the impurity models p9|. On the other hand, comparing the results 
for Ug = 8 and Ug = 10, there are no significant changes as long as z^ > 0.01 |^. This 
means (see Fig. ||) that the overall layer and U dependence is predicted reliably. We 
also believe that the finding of a common critical interaction is not an artifact of the ED 
approach since this is made plausible by the linearized DMFT. 

At the critical interaction the metallic solution continuously coalesces with the insu- 
lating solution that is found for U > f/c,buik- The insulating solution persists down to 
another (common) critical interaction strength f/c,i < f4,buik (we find t/c,i ~ 11.5). In the 
coexistence region, however, it is thermodynamically irrelevant. For details we refer to 
Refs. i, H, [1|, H . 



7.1 Modified intra- layer surface hopping 

A modification of the model parameters at the very surface may strongly affect the quasi- 
particle weight. As in Sec. ^ we first consider a modified hopping within the top layer: 
til ^ t. 

Fig. 1^ gives an overview for fixed Coulomb interaction U = 10. The above-mentioned 
tendency towards an insulating surface is enhanced when tn is decreased. The top-layer 
quasi-particle weight quickly decreases but even for tu = it does not vanish completely. 
For til > t one can see the opposite trend. Za=i increases with increasing tu. In the limit 
til I— * C)0 it approaches its non-interacting value Za=i = 1. For tu = lOt the low-energy 
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Fig. ^: The same as Fig. ^ but on an enlarged scale. 



electronic structure is almost perfectly decoupled. In the top surface layer there is a 
quasi uncorrelated motion of the electrons {za=i = 0.98 at U/tu = 1). The rest system, 
however, remains to be a strongly correlated Fermi liquid. 

For the sub-surface layers, the dependence of the quasi-particle weight on tu is com- 
paratively weak. Fig. |^ shows z^ for a > 2. On the enlarged scale in Fig. ^ there is still 
a considerable tu dependence of Za=2 (second layer). For a t— > 00, however, i. e. with 
increasing distance to the surface, this dependence diminishes: The bulk quasi-particle 
weight obviously cannot be affected by the surface modification of the hopping parameter. 
We also notice that there is a nearly constant quasi-particle weight for tu ~ 0.8t and all 
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a > 3. 

For fixed tn one finds an oscillating layer dependence of Za- This is demonstrated in 
Figs. H and |lO]for tu = t and tu = lOt. For the strongly perturbed system with tu = 10, 
the layer dependence is somehow irregular in the surface-near region, oscillations do not 
build up until a > 5. In both cases the oscillation is strongly damped. For a = 13 we 
have Az/z ~ 2 ■ 10~^. Thus, for a film with thickness d = 25, the quasi-particle weight 
is nearly constant at the film center. Furthermore, the differences between the uniform 
and the perturbed system become smaller and smaller with increasing distance to the 
surface. The observed oscillations can be traced back to oscillations of the free layer DOS 
at the Fermi energy. It is well known |^ that the presence of the surface gives rise to a 
layer-by-layer oscillation of pa{E = 0) for U = 0. For the present case (local self-energy, 
manifest particle-hole symmetry, metallic phase), the density of states at the Fermi edge 
is unrenormalized by the interaction (see Eq. (|TOD and Ref. |23[). The same oscillation 
is thus found for pa{E = 0) at any U < f/c,buik and will also lead to oscillations of the 
low-energy part of the Green function and thereby to oscillations of the low-energy part 
of the self-energy. Finally, the oscillating behavior of Za for U = 10 shows that we are well 
below the critical point: For U close to f/c,buik we expect a monotonous behavior from the 
hnearized DMFT (see Fig. |). 

Let us now tackle the question of surface phases. The scenario of an insulating surface 
coexisting with a metallic bulk was excluded by the linearized DMFT. The same is found 



by the numerical evaluation of the DMFT: Fig. 11 shows the layer-dependent quasi- 
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Fig. |ll|: U dependence of Za when the intra-layer surface hopping is 
switched off. 



particle weight for tn = where the strongest suppression of Za=i is expected. In fact, 
the top-layer quasi-particle weight quickly decreases as a function of U and, compared 
with the bulk value, becomes very small above U ^ 6. However, we find a non-zero weight 
in the top layer up to U = f/c,buik which implies Uc^surf = f^c,buik- Between U ^ 6 and 
U = Uc^hnik we may speak of an induced metallic surface according to the discussion in 
Sec. I ' 

The linearized DMFT predicted a metallic surface on top of a Mott-insulating bulk 
to be possible for tu > tJ5/A. We choose tn = 1.5t for the numerical calculation to be 
well above this threshold. Fig. ^ proves that two different critical interactions are found 
indeed. Over the whole U range considered, the top-layer quasi-particle weight is strongly 
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enhanced compared with the bulk and is finite also at f/ = ?7c,buik where the bulk weight 
vanishes. Note that Za=\{U^ is continuous at [/ = t/c,buik- The top-layer quasi-particle 
weight approaches zero at [/ = t/c,surf = 20.0 which marks the surface transition point 
while the extraordinary transition takes place at t/ = f/c,buik = 16.0. 

Evaluating the analytical formula for the surface critical interaction (^81) for the present 
case, we get f/c,surf = 18.2 which agrees well with the numerical result if one takes into 
account that also for the bulk critical interaction the linearized theory yields a somewhat 
smaller value. We also expect that t/c,surf (as ^c,buik) is overestimated by the ED due 
to finite-size effects p3|, |^. While finite-size effects prevent a precise determination of 
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the critical interactions, they are irrelevant concerning the very existence of the metallic 
surface phase. Even for U well above t/c,buik, the top-layer quasi-particle weight is still 
larger than Za=i = 0.01, and thus the ED for n^ = 8 is still able to resolve the energy 
scale set by the width of the Kondo-type resonance at the surface. 

Since the low-energy surface excitations cannot propagate into the bulk for U > t/c,buik 
but are reflected at the bulk Hubbard gap, the Kondo resonance represents a true surface 
state. Therefore, its amplitude must decay exponentially with increasing distance to the 
surface. Fig. O shows that some weight is induced in the sub-surface layers which rapidly 
decreases. 

The surface transition is also found as a function of tn for fixed interaction U > t/c,buik- 
Fig. |13] shows the numerical results for U = 18. When tn > tn c = 1.33t the surface 
becomes metallic. The critical value may be compared with tu^c = l-48t which is obtained 
by solving Eqs. (|T6|) and (|23| ) for tn. 



For the tn range considered in Fig. |T3|, the number of layers in the slab d that is 
necessary to simulate the actual surface can be lowered down to d ^ 5: We performed 
calculations for different thicknesses d; there are hardly any differences between the results 
for Za at the surface as long as d > 5. This is interpreted as follows: Since the coherent 
part of the bulk spectrum has disappeared for U > Uc,hu\k, the surface electronic structure 
is essentially decoupled from the bulk in the low-energy regime. The decoupling at the 



low-energy scale is indicated by the rapid decrease of Za with increasing a (see Fig. |T3D. 
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Contrary, on the high-energy scale set by the charge excitations, bulk and surface modes 
cannot decouple. There is always a finite energetic overlap of the bulk and the surface 
DOS since tn mainly changes the effective widths but not the positions of the Hubbard 
peaks in the surface DOS. The effect of the Hubbard bands on the low-energy features, 
however, seems to be rather weak since otherwise a change of d would lead to significant 
changes in the surface low-energy electronic structure by indirect coupling between low- 
and high-energy surface excitations and high-energy surface and bulk excitations. The 
surface Kondo resonance in the metallic surface phase is spatially confined to the first few 
layers and energetically isolated from the surface Hubbard bands. 



7.2 Modified inter-layer surface hopping 

A complete decoupling between the top layer and the rest system is obtained for vanishing 
inter-layer surface hopping ti2 = 0. Fig. |14| shows the layer-dependent quasi-particle 
weight as a function of ti2. While for f/ = 10 we noticed an oscillating layer dependence 
for uniform hopping (Figs. ^ and |10D, there is a monotonous layer dependence for U = 12 
(Fig. [1^, for ti2 = 1). This is the typical behavior when the system is close to criticality 
as has been noted before (cf. Ref. p9[ and the discussion of the analytical z^ profiles in 
Sec. ^. The layer dependence remains to be monotonous for ti2 ^—>- 0. For ti2 = we 
have essentially two independent systems. The isolated top layer is still metallic. In the 
rest system the a = 2 layer represents the new top layer, the a = 3 layer becomes the 
first sub-surface layer, and so on. This implies that the value of Za for tu = must be 
equal to the value of Za-i for ti2 = 1. These relations are indicated by the dashed lines 




Fig. |l^: The same as Fig. |l^ but for ti2 > t = 1. Inset: asymptotic 
behavior of Za- 
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in Fig. |T^. They represent a non-trivial check of the numerics. 

An effective separation into subsystems is also observed in the opposite limit of a 
strongly enhanced inter-layer surface hopping. Fig. [l^ shows that Za=i and Za=2 approach 
their non-interacting values while for a > 3 the quasi-particle weight changes only slightly 
as ti2 ^—>- cxD. In the low-energy regime the electronic structure of the first two layers 
decouples from the rest system. The value of Za for all a > 3 approaches the value of 
Za-2 for ti2 = 1 (see inset). 

A somewhat artificial realization of an insulating surface phase on top of a metallic bulk 
can be obtained for tu i— *> by choosing Uc,d=2 <U < t/c,buik, where Uc,d=2 is the critical 
interaction of the two-dimensional layer. For U > Uc,d=2 the top layer must become 
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insulating when it is decoupled from the rest system {tu = 0). This is demonstrated in 
Fig. 0. The figure also shows that the top layer becomes metallic (with a very small 
quasi-particle weight) as soon as an arbitrarily small inter-layer hopping is switched on. 
Finally, Fig. p!7| shows the extraordinary and the surface transition for fixed ti2 = 3.0. 
The surface critical interaction can be read off to be Uc^surt = 23.8 while the linearized 
DMFT with ?7c,surf = 21.7 [Eq. (0)] again predicts a slightly smaller critical value. 

7.3 Modified surface Coulomb interaction 



We finally discuss the modification of the surface Coulomb interaction Ui. Fig. |TB| shows 
the quasi-particle weight Za for a = 1,2,3 and 2;buik as a function of Ui/U where U is 
fixed at U = 10. On decreasing Ui {Ui < U), zi quickly increases, and for [/i h-> it 
approaches the non- interacting value zi = 1. 

For enhanced Ui > U the top-layer quasi-particle weight is decreased but remains to 
be finite even for large large values of f/i, i. e. we again find an induced metallic surface. 
Asymptotically, however, the top-layer weight approaches zero: zi i— > for f/i \—>- cxo. 
In this limit the low-energy resonance in the top-layer DOS essentially disappears and 
a large Hubbard gap ~ f/i opens. This implies that - in the low-energy regime - the 
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Fig. |T9|: U dependence of Za for enhanced surface Coulomb interaction 
Ui/U = const. 



sub-surface (a = 2) DOS for Ui ^-^ oo must become identical with the Ui = U surface 
(a = 1) DOS. For f/i i— > oo we thus expect ^2(^1 *-^ 00) = Zi{Ui = U) and consequently 
ZaiJJi I— > 00) = Za-iiJJi = U) for all a. In fact, this "shift" of the surface by one layer 
can be seen in Fig. |18| and in the inset: For Ui/U = 5 only small differences still remain 
between Za and Za-i{Ui = U). 

The "shift" a ^^ a — 1 also implies that the oscillating layer dependence of z^ for 
Ui = U must be reversed for Ui ^-^ 00: Minima are replaced by maxima and vice versa. 
This partly explains that between (at Ui/U ~ 1.2) the quasi-particle weight is nearly 
layer independent. 

According to the linearized DMFT, a metallic surface on top of an insulating bulk can 
be found if Ui < Qt^Jl [Eq. (H)] and U > Uc,hnik = StVe. If we fix the ratio Ui/U = 0.5 
and vary U, the surface transition should occur at t/i,c,surf = 12.1 (set U = 2Ui^c,sut{ in Eq. 
(|3^ ) and solve for t/i,c,surf)- The result of the numerical solution of the DMFT equations 
is shown in Fig. ^. Again, the numerically obtained value, t/i,c,surf = 0.5 • 26.3 = 13.15, 
is somewhat larger than the analytical prediction - the discussion is the same as for the 
modified surface hopping. 

At the extraordinary transition U = f/c,buik; the quasi-particle weight in the top- (not 
shown in Fig. |l^) and in the first sub-surface layers are smooth functions of U. This is 
contrary to the results found within the slave-boson theory [^ where the band-narrowing 
factor in the top layer shows up a discontinuous derivative at f/ = f/c,buik- We believe, 
however, that this is an artifact due to incorrect boundary conditions. Namely, only the 
first two surface layers are treated as "free" in the self-consistent calculation while already 
the a = 3 layer is assumed to be bulk-like in Ref. [^. As is known from the mean- field 



theory of localized spin models 



such boundary conditions may result in artificial 



singularities at the extraordinary transition. 
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8 Conclusion 

We have investigated the (Mott) metal-insulator transition at surfaces within the frame- 
work of the semi-infinite Hubbard model at half- filling and T = 0. Basically, two approx- 
imations have been used: 

First, the self-energy functional has been assumed to be reasonably local. This approx- 
imation sets the basis for the dynamical mean-field theory: The semi-infinite Hubbard 
model is self-consistently mapped onto a set of indirectly coupled impurity models corre- 
sponding to the inequivalent layers parallel to the surface. With the usual scaling of the 
(intra-layer and inter-layer) hopping, the approach becomes exact in the limit of infinite 
spatial dimensions. It has been shown that there are non-trivial surface effects even for 
D = oo. Mainly, however, the DMFT has been used as a (mean- field) approach to study 
the D = 3 low-index surfaces of the simple-cubic lattice. 

Second, for the approximate solution of the impurity models, we have used the exact 
diagonalization of finite systems. The ED method allows to deal systematically with a 
large number of geometries and model parameters. However, the method cannot access 
the very critical regime for the Mott transition because of errors due to finite-size ef- 
fects. Directly at the critical point we have alternatively considered a simplification of 
the mean-field equations (linearized DMFT). This analytical approach is also approxi- 
mate. However, a convincing qualitative and (as far as can be judged) also quantitative 
agreement with the numerical ED results has been found. 

Referring to the points mentioned in the introduction, our results can be summarized 
as follows: 

1. The metal-insulator transition in the bulk of the semi-infinite system occurs exactly 
at the same critical interaction t/c,buik as for the infinitely extended system: t/c,buik = 

2. There is a non-trivial layer dependence of the quasi-particle weight, even (asymptot- 
ically) at the critical point. The Za profile strongly depends on the model parameters 
at the surface, e. g. the hopping within the top layer tn, the hopping between the 
top and the sub-surface layer tu and the top-layer Coulomb interaction Ui. There 
is a qualitative change of the profile if certain critical values, tn^c, ti2,c and [/i,c, are 
exceeded. These critical values are found to be of a realistic order of magnitude. 

3. For uniform model parameters the top-layer quasi-particle weight Zi is smaller than 
the bulk value -Zbuik since a reduced surface coordination number implies correlation 
effects to be effectively stronger at the surface. For interactions well below Uc-, there 
is always an oscillating layer dependence of the the quasi-particle weight. With 
increasing distance to the surface (a i-^ oo), this oscillation is strongly damped. 
In the critical regime, on the contrary, Za monotonously increases with increasing 
a, and finally, for U = Uc the critical profile is linear: Za oc a. For uniform model 
parameters there is a finite weight in the top layer {zi > 0) for U < Uc only, i. e. only 
when the bulk is metallic. The transition at Uc is termed the "ordinary transition" . 

4. For a sufficiently strong modification of the surface model parameters (tn > tn^c, 
ti2 > ti2c, Ui < Uic), the surface becomes metallic below a critical interaction 
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f^c,surf > f^c.buik ("surface transition"). For f/c,buik < U < f/c,surf, the quasi-particle 
weight exponentially decays from its maximum value zi at the surface towards zero 
in the bulk. At U = f/c,buik the bulk undergoes the transition to the metallic 
state ("extraordinary transition"). The top-layer quasi-particle weight is a smooth 
function of U even at f/ = f/c,buik- 

5. The transition at f/ = f/c,surf = f^c.buik is called the "special transition". Here the 
critical profile of the quasi-particle weight is flat Za = -Zbuik = const, (at least for 
a > 2). In this situation the effect of missing neighbors at the surface is compensated 
by the change of the surface model parameters. 

6. There are two critical exponents that are merely related to the critical interactions, 
the "shift exponent" A^ and the "crossover exponent" 0. They describe the trend 
of the Uc for a film of finite thickness d in the limit d \^ oo and the trend of the 
surface critical interaction for the semi-infinite system near the special transition, 
respectively. Within the linearized DMFT one finds A^ = 2 and = 1/2. 

7. For any realistic choice of the model parameters, a metallic bulk induces a metallic 
surface with zi > 0. Thus, a Mott-insulating surface of a correlated metal is impos- 
sible. There are essentially two more or less trivial exceptions: The first is the static 
decoupling of the top layer for t^ = at an interaction strength that is smaller 
than f/c,buik but larger than the critical interaction of the two-dimensional system. 
The second is a dynamical decoupling which occurs for infinite surface interaction 
Ui \-^ CO. Here the top-layer quasi-particle weight vanishes asymptotically. 
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